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, Abstract The negative effective magnetic pressure instability operates on scales 

C^ ' encompassing many turbulent eddies and is here discussed in connection with 

the formation of active regions near the surface layers of the Sun. This instability 

is related to the negative contribution of turbulence to the mean magnetic 

^«0 I pressure that causes the formation of large-scale magnetic structures. For an 

isothermal layer, direct numerical simulations and mean-field simulations of this 

^^ I phenomenon are shown to agree in many details in that their onset occurs at 

r^ • the same depth. This depth increases with increasing field strength, such that 

the maximum growth rate of this instability is independent of the field strength, 

provided the magnetic structures arc fully contained within the domain. A linear 

stability analysis is shown to support this finding. The instability also leads to a 

O . redistribution of turbulent intensity and gas pressure that could provide direct 

*-J ' observational signatures. 
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Active region formation in the Sun is traditionally thought to be a deeply 
rooted phenomenon, because their size (^ 100 Mm) is much larger than the 
naturally oc curring scale s in th e surface layers of the convection zone {'^ 1- 
10 Mm); see iGolub et al.l (|l981l ). They are also long-lived (~ 1/2 year), which 



seems unnaturally long if associated with the near-surface layers (40 Mm depth) 
where typical time scales are about a day. On the other hand, a deeply rooted 
formation scenario for active regions has the problem that the azimuthal pattern 
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speed of active regions does not match the angular velocity at great depth. Other 
difficulties concern the strong field strength inferred for the tachocline to explain 
the observed tilt angles and the fact that magnetic structures expand tremen- 
dously during their ascent. These and several other argu ments have led to the 



consideration of solar activity as a shallow phenomenon; see lBrandenburd (|2005f ) 
for details. As a possible mechanism for producing magnetic flux concentrations 
of the form of active reg ions, the negative effective magnetic pressure insta- 



bility has been discussed (IKleeorin et al 



1994; Klecorin . Mqnd. an d Rogachevskii , 



198 



199 



1990t Kleeorin fc Rogachevskn . 



Rogachevskii fc Kleeorinl . 12007 : 



Brandenburg et al.l . 12010 . 2011 ). Of course a magnetic field B always gives 



rise to a positive magnetic pressure, B'^/2fio, where /io is the vacuum per- 
meability. In a turbulent medium, however, magnetic fields also suppress the 
turbulence and thus decrease the turbulent pressure pu^/3, and modify the 
pressure caused by magnetic fluctuations b^/6/^o- Here, u and b are velocity and 
magnetic fluctuations, p is the density, ^o is the vacuum permeability, and the 
coefficients in the turbulent ffuid and magnetic pressure are given for isotropic 
turbulence. Magnetic fluctuations can be due to both small-scale dynamo action 
as well as tangling of a large-scale field, B. The total field is thus B = B + b. 
The sum of both effects, pturb = pu-^/i + {b'^)/6po, is positive definite, but 
it depends on B, and pturb tends to decline as i? = \B\ increases. Indeed, 
Pturb = '2Ek/3 — (6^)/6/^o, where Ek = pu'^/2 + (b^)/2/io ~ const, so that 
the change of the turbulent pressure is negative (Spturh < 0) when the magnetic 
fluctuations are generated by tangling of the mean magnet ic field by the v e locity 
fluctuations at the expan se of turbulent kinetic energy (IKleeorin et al.l . Il990l : 
Brandenburg et al.l . 120111) . Thus, we write 



Pturb(B) = Pturb(O) - gp(S)B'/2Mo, (1) 

where pturb(O) is the turbulent pressure at zero mean field. The pressure pturb(O) 
only includes those contributions from b^ that are associated with small-scale 
dynamo action, but not those resulting from the mean field. The relevant mag- 
netic pressure in the evolution equation for the mean flow, U, is then not just 
B /2fiQ, but it is affected by the B dependence of pturb, i-e., it depends on 

Pturb(S) + B^/2^Io = Pturb(O) + [1 - (7p(S)]s'/2Aio, (2) 

which is also still positive, but 1 — q-p{B) may well become negative, which is 
what we call a negative effective magnetic pressure. Consequently, the expression 

Pcff = (1 - gp)BV2Aio (3) 

is referred to as the effective magnetic pressure. In addition, there is also the 
gas pressure Pgas- Once the effective magnetic pressure drives a mean flow, the 
gas density changes, and as a consequence the gas pressure, so as to re-establish 
approximate total pressure balance. Therefore, pgas and p will also depend on 

In the presence of gravity, the properties of magnetic buoyancy are drastically 
altered by a negative effective magnetic pressure. In the following we illustrate 



SOLA: paper.tex; 25 February 2013; 10:19; p. 2 



Negative effective magnetic pressure instability 



how this can lead to an instabiUty. Since the flow velocities are highly subsonic, 
we can make the anelastic approximation for low Mach numbers, i.e., V-pC/ = 0. 
This leads to V • 17 + C7 • V Inp = 0, or 

V-C7=^, (4) 

tip 

where we have used the density scale height Hp, so that Vhip = (0,0, —l/Hp). 
This equation shows that a downward motion C/^ < leads to a compression, V • 
[/ < 0. This enhances an applied field locally. We consider an applied equilibrium 
magnetic field of the form (0, Bq, 0) and the mean field has only a y component, 
i.e., B = (0, By{x, 2;), 0), so we have 






-ByV ■ U, (5) 



where D/Dt = d/d+ C/ • V is the advective derivative. Note that for a magnetic 
field with only a y component, but d/dy = 0, there is no stretching term, so there 
is no term of the form B ■ VC/. Using Equation (j4|), and linearizing Equation ([5]) 
around U = and B = Bq, we have 

where subscripts 1 denote linearized quantities. The vertical velocity perturba- 
tion Uiz is caused by magnetic buoyancy. Assuming total pressure equilibrium, 
Pgas + PcS = const, we see that an increase in the effective magnetic pressure 
causes a decrease is the gas pressure, i.e., (5pgas = —5pcS, just like in the regular 
magnetic buoyancy instability. Therefore, the Archimcdian buoyancy force is 

2 

5p 6pg^s Spctf dpcs SB 

P Pgas Pci dB P^i 

where we have used Pgas ~ pc^ for an isothermal gas. In the regular magnetic 
buoyancy instability, without turbulence effects, we have 2/iodpoff/di3'^ = 1. 
In the domain where the negative effective magnetic pressure effect causes 
dpeff/dB^ to be negative, a magnetic field enhancement leads to a further 
reduction of the local pressure, which is compensated by horizontal inflows, 
increasing density (and field strength), making this fluid parcel heavier, causing 
it to sink. Inversely, a local field reduction causes outfiows and rises until it 
reaches the region where this feedback reverses. Thus, the instability loop is 
closed by considering the momentum equation in its linearized form 

dUiz dpcff 2BoBiy 

—5^ = ^=2" 2 9- 8) 

dt dB PCs 
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Using Cg/g = Hp for an isothermal atmosphere, we then find the dispersion 
relation for the growth rate A of the resulting instability 



-t 



-J~2iiop',s-Vtk\ (9) 

where v^ = Bo/y/JiQp is the Alfven speed and 

2 2 

2/io Pcff = 2/^0 dpcff /dS = 1 - gp - dqp/d In _B (10) 

is twice the derivative of the effective magnetic pressure. We have here also 
included the effects of turbulent magnetic diffusivity rjt and turbulent magnetic 
viscosity ft, assuming vt/vt = 1- Here, k is the effective wavenumbcr. A proper 
derivation of the growth rate of the instability, but again without including 
turbulent magnetic diffusivity and turbulent magnetic viscosity, is given in 
Appendix [X] 

The negative contribution of turbulence to the mean magnet ic pressure and 
the resulting large-scale instability has been predicted long ago (IKleeorin et al 



IQQOtlKleeorin. Mond. and RogachevskiiLllQQalKleeorin fc Rogachevskiilll994l) 



However, this ins t ability has been , detec ted in DNS only recently 
( Brandenburg et all I2OIII : iKemel et al.l . l2012al) . This large-scale instability is 



called the negative effective magnetic pressure instability, or NEMPI, for short. 

Equation (jH]) demonstrates that stronger stratification and thus a smaller scale 

height leads t o an increased growth rate of the instability. This was qualitatively 

confirmed by iKemel et al.l ( 2012bl ) . Using numerical solutions of the full mean- 



field equations, they found furthermore that the maximum growth rate of the 
instability is actually independent of va- This seems to be at odds with the 
Equation ^. To understand this, we use the following fit formula for gpi 



PI 



'^P(^) - ^^^' (11) 



where /? = B / B^^ and B^q = y/J^opu^-ms is the equipartition field strength. Thus, 
for Pi, ^ /3 :^ I3p, we have 

A«/3.^-r;,fc^ (12) 

Bp 

so the growth rate is indeed independent of the imposed field strength. 

In a mean-field model, Wmis is normally expressed in terms of rjt = Wmis/Sfef , 
so Equation p^ turns into 

' -3/3.^-1, (13) 



rjtk'^ kHp 

which illustrates immediately the importance of large enough scale separation, 
i.e., large enough values oi kf/k. 

The purpose of this paper is to show that NEMPI can work over a range 
of different field strength. S uch a result was rece ntly predicted using the mean- 
field simulations (MPS) bv iKemel et al.l ( 2012br ). We shall also investigate the 
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close connection between mean field and the resulting effective magnetic pres- 
sure in the plane perpendicular to the mean field. Here we focus on a series of 
simulations with different field strengths, but for a fixed value of the magnetic 
Reynolds number and fixed value of the scale separation ratio. For a study of 
the dependence on magnetic Reynolds number and on scale separation ra tio, but 
fixed field strength, we refer to the recent work of iKemel et al.l (|2012a[ ). In the 



following, we discuss first direct numerical simulations (DNS) of NEMPI and 
turn then to proper mean-field simulations (MFS). We begin with a simplistic 
illustration of the nature of NEMPI. 



2. Vertical profile of efftective magnetic pressure 

The first successful DNS of NEMPI have been possible under the assump- 
tion of an isothermally st ratified layer with an isothermal equation of state 



(jBrandenburg et al.l 120111) . Much of the same physics is also possible in adi- 
abatical ly stratified layers, but NEMPI was found in th is case only in mean-field 
models ( Brandenburg et al.l . l201CltlKapvla et al.l . 120121 ). The isothermal case has 



conceptual advantages that help us understanding better the underlying physics 
of this instability. We make use of this advantage in the present paper, too. 

In most of the isothermal setups studied so far, the rms velocity is only weakly 
dependent on height, so the z variation of i?cq was only caused by that of p = 
pQ exp{— z / H p) . This allows us then to plot the effective magnetic pressure. In 
the following, we introduce the quantity 

r,s{p) = \\i-q^mP\ (14) 

which is the effective magnetic pressure normalized by the local equipartition 
field strength Scq, i.e., V^s = Pes/B'^^. Since /3 = /3{z) = /3o exp(z/27Jp) in- 
creases with z, Veg{z) is small at large depths, reaches a negative minimum at 
some depth, and then becomes positive and equal to /3^. In Figure [1] we show 
vertical profiles oiVeS, dT'off/d/S^, and /3(— 2d'Peff/d/3^)^/^ for the fit parameters 
gpo = 20 and /3p = 0.167 derived later in this paper, and the three field strengths 
/3q = Bo/BcqQ = 0.05, 0.1, and 0.2 within the z range from ~~n to tt, which 
is also consistent with the DNS and some of the MFS discussed below. Here, 

-SoqO = Soq(O). 

Notice first of all that all three curves of T^cff have minima with left fianks 
(negative slopes) within the domain. As the imposed field is increased, these 
curves shift downward (smaller values of z). Thus, we should expect the peak of 
the instable eigenmode to appear somewhere along the left flanks of these curves 
and that these peaks move further down as the imposed fleld is increased. This 
is qualitatively reproduced by the DNS and MFS discussed below, except that 
the location is consistently a certain distance below the position where the left 
fianks have their steepest gradient. On the other hand, as is evident from the 
middle panel of Figure \T\ the largest value of dPeff/d/?^ is always achieved at 
the bottom of the domain. However, the growth rate of NEMPI has still a factor 
proportional to va — u^msP in front of it; sec Ec^uation ^. This then confines the 
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Figure 1. Profiles of -Peff, AV^slAfi^-, and /3(-2d'Peff/d/3^)^/^ for fit parameters Qpo = 20 
and /3p = 0.167, and the three field strengths _Bo/BeqO = 0.05, 0.1, and 0.2 within the z range 
from — TT to TT, which is consistent with some of the models discussed below. The vertical lines 
of similar line types give the location where the unstable cigenmodes reaches its peak. 



instability to a narrow strip within the domain. In the third panel of Figure[T]we 
plot therefore also /3(— 2d7'off/d/3^)^/-^, and their extrema are now only slightly 
above the location where DNS and MFS show a peak in the eigenfunction. The 
reason for the remaining discrepancy is not well understood at present. 
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3. Onset and saturation of NEMPI in DNS 

3.1. Isothermal setup in DNS 



Following the earlier work of iBrandenburg et al.l (1201 if ) and iKemel et al. 



(|2012al ) , we solve the equations for the velocity U, the magnetic vector potential 



j4, and the density p, 



P^ = -^s^P + Jy-B + p{f + g) + V- (2i.pS), (15) 



Df 



dA 



^UxB + jyV^A, (16) 

I ^ -V.pt/, (17) 

where i' is the kinematic viscosity, 77 is the magnetic diffusivity due to Spitzer 
conductivity of the plasma, B = Bq+'VxAis the magnetic field, Bq = (0, Bq, 0) 
is the imposed uniform field, J = V x B/fiQ is the current density, /zq is the 
vacuum permeability, S^ = ■^(Ui,j + Uj,i) — h^ij^ ■ U is the traceless rate of 
strain tensor, and commas denote partial differentiation. The forcing function 
/ consists of random, whitc-in-timc, plane non-polarized waves with a certain 
average wavenumber. The turbulent rms velocity is approximately independent 
of z with Urms — (u,^)^/^ w 0.1 Cg. Thc gravitational acceleration g ~ (0, 0,— g) 
is chosen such that kiHp = 1, so the density contrast between bottom and top is 
exp(27r) sa 535. Here, Hp = c^/g is the density scale height. We consider a domain 
of size Lx X Ly x L^ in Cartesian coordinates {x^y,z), with periodic boundary 
conditions in the x and y directions and stress-free perfectly conducting bound- 
aries at top and bottom {z = ±Lz/2). In all cases, we use a scale separation 
ratio ki/ki of 30, a fluid Reynolds number Re = Umis/^^^f of 18, and a magnetic 
Prandtl number Pr^/ ~ v jr] of 0.5. In our units, [Lq = \ and Cg = 1. The value of 
-Bq is specified in units of the volume averaged value, Boqo = y/fJ-oPo Urms, where 
Po — (p) is the volume-averaged density, which is constant in time. In addition 
to visualizations of the actual magnetic field, we also monitor By, which is an 
average over y and a certain time interval At. Since the simulations arc periodic 
in the x and y directions, we sometimes shift the images such that the peak field 
strength of NEMPI appears in the middle of the frame. 

The simulations are performed with thc Pencil Code|^ which uses sixth- 
order explicit finite differences in space and a third-order accurate time stepping 
method. We use a numerical resolution of 256^ mesh points. 

3.2. Resuhs 

In FigurelHwe demonstrate that NEMPI can work over a range of field strengths. 
As we increase thc strength of the imposed field, NEMPI develops at progres- 
sively greater depth. This result was recently obtained for MPS, but is now 



^http : //pencil-code ■ googlecode ■ corn] 
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Figure 2. By from DNS for three values of the imposed field strength at the end of the linear 
growth phase of NEMPI for Ksm = 18 and Prjvf = 0.5. 
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Figure 3. Growth of the large-scale field strength _Bi at the center of the magnetic structure 
for three field strengths. 



for the first time demonstrated in DNS. Figure [3] shows that the growth of the 
large-scale field Bi of the magnetic structure is similar for three different field 
strengths. Here, Bi has been determined by taking the maximum value of the 
mean field in the neighborhood of the position where the flux concentration 
later develops. Note that there is a range over which Bi grows approximately 
exponentially, independent of the value of Bq. 

In FigurelHwe show By at early, intermediate, and late stages of the saturation 
process (left), and compare with visualizations of VeS at the same times. Here, 
'Pcff = ■^(1 ^ Qp)P'^i where (7p(/3) with /? = B/Bcq is evaluated from 



gp = -2AU^JB\ 



for JBo = (0,Bo,0), and 



AH, 



-p{u^^ul)+'^{b^-bl)-{b^-bl), 



(18) 



(19) 
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Figure 4. By and VcS from DNS at tlirce times sliowing the descent of the potato sack 
feature for Rej\/ = 18 and Pr^i/ = 0.5. 



is applied to the xx component of the total stress from the fluctuating velocity 
and magnetic fields. In Equation (|19p no summation over the repeated index i 
is assumed. 

In Figure |31 blue shades correspond to low values of VcS and occur around 
the minimum line (marked in white) where Ves = T'min- As time progresses, low 
values of PeS are also found at greater depth as the magnetic flux concentration 
descends. The fact that there is a clear spatial correlation between By and T^off 
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iC 1 X 



Figure 5. Vectors of U together with a color/grey scale representation of By from DNS at 
a late time for Rej\/ = 18 and Prj\/ = 0.5. 



provides strong evidence that the interpretation of the formation of structures 
in the stratified turbulence simulations in terms of NEMPI is indeed the correct 



one. 



The descendi ng structures have previ ously been referred to as potato 
sack structures (JBrandenburg et al.l . l201l[) . because of their widening cross- 
section with greater dept h. When such structures were first seen in MFS 
( Brandenburg et al.l |2010|), they were originally thought to be artifacts of the 
model that one would not ex pect to see in the Sun. H owever, such structures 
were later also found in DNS (JBrandenburg et al.l . 120111) , highlighting therefore 
the strong predictive power of MFS. 

Visualizations of the resulting mean flow U are shown in Figure [S] as 
vectors. The flow shows a convergent shape toward the magnetic struc- 
tures. It is interesting to note that such convergent flow structures are 
now also seen in local heliosei s mic f low measurements around active regions 
( Hindman. Haber. and Toomrel 120091 ) . In this connection it is instructive to 
discuss the somewhat peculiar shape of such a structure that widens as it 
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descends. Normally, in a strongly stratified atmosphere, descending structures 
get compressed and become nar rower, but this is not seen in the present visu- 
alizations. As already argued in [Brandenburg et al.l ( 2012 ). this is because the 
boundaries of these structures do not coincide with material lines, so the mass 
is not conserved inside them and can leak through the boundaries. Indeed, these 
structures gr ow as they descend, and may becorn e amenable to helioseismic 
detection; cf. Illonidis. Zhao, and KosovichevI (J201lh . This phenomenon is well 
known in the d escription of turbulent p lumes as a model of turbulent downdrafts 
in convection (jRieutord fc Zahnl . 119951 ) . Such structures are known to widen as 
a result of entrainment. The sinking behavior of these apparently disconnected 
flow structures can be explained as follows: while inflows dominate downflows 
throughout the whole lifespan of the field concentration, in the initial stage the 
former can drag in a large fraction of the surrounding magnetic field, overcom- 
pensating the losses by downflows. However as the environment gets depleted, 
this dynamical balance shifts and the structures start moving downwards. 



3.3. Mean-field coefficients from DNS 



In earlier work bv lKemel et al.l (|2012b[ ). the parameters gpo = 40 and /3p = 0.05, 
corresponding to /3^ = 0.3 2, w ere used. Those value s are compatible with work by 
Brandenburg et al.l ( 2012 ) and lKemel et al.l ( 2012ai ). However, in the present case 
we have a larger scale separation ratio, kf/ki ~ 30, for which these parameters 
have not yet been determined. In the Figure [5] we show the functional form 
of Vea{fi) for the present case with fcf/fci = 3 0, ReM = 18, and Pr/yf = 0.5. 



Here we have followed the method described in iBrandenburg et al.l (2012); see 
their Eq. (17). For the present model we find as fit parameters ^po = 20 and 
/3p = 0.167, which corresponds to /3* = 0.75. 



4. Comparison with MFS 



Recently, many aspects of NEMPI seen in the DNS have also been detected in 
MFS. Establishing the usefulness and limitations of MFS is important, because 
such models are easier to solve and allow one to explore parameters in regimes 
where DNS are harder to apply or have not yet been applied in the limited time 
span since the close correspondence between DNS and MFS was first noted. 

In the following we consider two-dimensiona l mean-field models , in which 
the presence of qs has no effect on the solutions (JKemel et al.l . l2012bl ). Further- 
more, we ignore other effects connected with the ani sotropy of the turbulence . 



These effects have p reviously been found to be weak (jBrandenburg et al, 
Kapvlaet all 120121) . 



2012; 



4.1. Isothermal setup in MFS 



In this section we solve the evolution equations for mean velocity [/, mean 
density p, and mean vector potential A, in the form 



QU 

— = -[/• VtZ-c^Vlnp 



■M 



^K, 



(20) 
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Figure 6. Peff {/?) for the DNS used in tliis paper with kf/ki = 30, Rci\/ = 18, and Pr^/ = 0.5. 






-U-Vp-pV- u, 



(21) 



dA 
'dt 



^UxB-{Tjt+v)J, 



where J^u is given by 

P^M = -iV[(l - qp)W] + B- V [(1 - qs)B] , 

and 

^K = {i^t + v) ("^^-U + i VV • U + 2SV Inp) 



(22) 



(23) 



(24) 



is the total (turbulent plus microscopic) viscous force. Here, S^- = ^{Uij + 
Uj^i) — hSij'V ■ U is the traceless rate of strain tensor of the mean flow. 

We approximate Qp by a simple profile that is only a function of the ratio 
/3 = \B\/Bcq. We use an algebraic fit of the form 



qpW) = - 



gpo 



/3V/?r 



(25) 



The function qp quantifies the impact of the mean magnetic field on the effective 
pressure force. 
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Bo/B.^o= 0-20 
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Figure 7. By from mean-field models for three values of the imposed field strength at the 
end of the linear growth phase of NEMPI. Here, gpo = 20 and /3p = 0.167, which corresponds 
to P^, = 0.75. 



Bq/B^o = 0.20 Bo/B,^a = 0.40 B^/B,^^= 0.80 

^ 1^ Wa 
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Figure 8. Similar to Figure [T] but for larger field strengths and vertical domain boundaries 
that a deeper down, so the magnetic field maxima of the instability fit better into the domain. 
Here, QpQ = 40 and /3p = 0.05, corresponding to /J* = 0.32. 



4.2. Aspects of the MFS 

We begin by showing By for three values of the imposed field strength at the 
end of the linear growth phase of NEMPI. The results are shown in Figures [7] 
and m for two different setups. In the former we use q^o = 20 and /?p = 0.167 
for the same z range (— tt < z/Hp < tt) as in the DNS, while in the latter 
we use Qpo = 40 and /3p = 0.05 for somewhat stronger fields an d a deeper z_ 



range ( < z/Hp < 27r), which is also the fiducial model used by iKemel et al. 
(|2012br ). In the former case the growth rate is w llH^/rjt while in the latter it 
is « 5.0H^p/T]t. 

Unlike the DNS, the MFS show that in the former series of models with qpo = 
20 and /?p = 0.167 the x extend is slightly larger than the optimal horizontal 
wavelength of the instability, because one sees that some of the structures begin 
to split into two (Figure [7]). This is not the case for the second model with 
qpO = 40 and /3p = 0.05_(Figure[8l). 

Next, we compare By with T'off = ^(1 ^ ^p)/?^- Again, there is a close 
correspondence between the By field and the resulting distribution of VqB] 
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Figure 9. By and Peff from mean-field models at three times. 



see Figure El Here, qp{f3) is evaluated from the assumed fit formula given by 
Equation (j25|) . Furthermore, there is a close correspondence between regions of 
enhanced magnetic field and enhanced density; see Figure [TUl 
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Figure 10. By compared with p for three different times. 



5. Conclusions 



The present work has demonstrated that NEMPI is able to concentrate the 
magnetic field into large patches encompassing the size of many turbulent 
eddies. The physics of this instability is a s t raightforwa r d extension of the 
usual magnetic buoyancy ins t abiliti e s (IParkeii 19661 1979t [Hughes fc Proctor . 
1988t .C attan eo and HugheJ . Il988l : IWissink et aiTbOOOt Tlsobe et all . 12005 : 
Kersale. Hughes, and Tobiasl . 120071 ) , except that the sign of the buoyancy force is 
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reversed in a regime of intermediate field strength. We have here re-examined the 
simple case of an isothermal layer in which NEMPI can in principle occur at any 
depth whose value is determined by the strength of the imposed field. Our new 
DNS have verified that the growth rate is indeed independent of the strength of 
the imposed field, provided the peak of the instability fits still comfortably within 
the domain. During the subsequent nonlinear evolution of the instability, the 
overall density stratification readjusts, allowing the magnetic field concentrations 
to move further down. It is important to reali ze that the resulting stru ctures are 
subject to significant turbulent entrainment (jRieutord fc Zahnl . 119951 ). so their 
boundaries are not closed. The agreement with corresponding mean field models 
is astounding and much more convincing than what has so been possible to 
demonstrate in mean-field dynamo theory. Mean-field models provide therefore 
a strong source of guidance when designing new setups for DNS. 

While NEMPI now begins to be fairly well understood for isothermal models, 
more work is required for non-isothermal ones. In that case, the density scale 
height is no longer constant and the degree of stratification is much stronger 
at the top than in deeper layers. The simple result that the instability can 
occur at any height, depending just on the strength of the imposed field, is then 
no longer valid. At the same time, there is another perhaps more important 
aspect, namely the possibility of other instabilities. One of them is connected 
with the suppressio n of turbulent convecti v e ener gy flux by the mean magnetic 
field. As shown by iKitchatinov &: Mazurl (2000|), this effect can also lead to 
magnetic flux concentrations and it may be sustained for much stronger magnetic 
field strengths, allowing thus the formation of structures in which the magnetic 
pressure becomes comparable to the ambient gas pressure. There may also be a 
connection with flux segre gation events seen in simulations of ma gnetoconvec- 
tion at large aspect ratios ( Tao et al.l . Il998t iKitiashvili et al.l . 120101 ) , which have 
already been shown to pr oduce bipolar regions in simulations with radiation 
transfer (jStein et aL . 2011 ). The s tudy of the possibility of producing sunspots 
similar to those of lRempell ( 2Qllal lbl). but without initial fiux structures, is now 
of high priority in the quest for solving the solar dynamo problem in terms of 
distributed dynamo models in which magnetic activity is explained as a surface 
phenomenon. 
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Appendix 

A. Growth rate of NEMPI 

In this Appendix we derive the growth rate of NEMPI neglecting for simplicity 
dissipation processes, using anelastic approximation for small Mach numbers 
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and assuming the density hight, Hp to be constant and jjlq 
the equation of motion in the following form 



dU{t,x,z) 

Ft 



Ptot 



Ptot 

pHp 



1. Let us rewrite 



(26) 



where z is the vertical unit vector, ptot = P + pcff is the total pressure (the sum 
of the mean gas pressure, p, and the effective magnetic pressure, Peff), and we 
took into account that mean magnetic field is independent on y, so that the 
mean magnetic tension vanishes. We also used an identity: 



Ptot 



Vpt, 



Ptot 



P pHp 

Taking twice curl of Equation (pS)) we obtain 

1 



i-(A-^)u.^ 



dt 



H, 



PH, 



A_LPcff, 



(27) 



(28) 



where A^ = A — V^ and we have used Equation Q. Introducing a new variable 
Vz = Vp Uz , we rewrite Equation (|28p for a new variable: 



^(A 
dt 



1 



AHl 



Vz 



1 ^ 

AxPcff- 
VpHp 



(29) 



Linearizing Equation (j29p . and using the linearized induction Equation ([5]) we 
arrive at the following equation: 



dt^ \ ml 



Vz{t,x,z) 



2131 dpcff 



-pHj d/32 



A^V^,, 



(30) 



Po 



where dpcft/d/^^ = 1(1 - gp - dqp/dln/?^)^^//?^ /? = B/B,^{z) and /3o = 
Bo/Bcqo- It follows from Equation ((30)) that a necessary condition for the large- 
scale instability is 



dp, 



'eff 



d/32 



<0. 



(31) 



Po 



For instance, in WKB approximation when k^ Hp 3> 1, i.e. when the charac- 
teristic scale of the spatial variation of the perturbations of the magnetic and 
velocity fields are much smaller than the density hight length Hp, the growth 
rate of the instability reads 



PokA 



fp^Hpk 



-2 



dpcff 



d^2 



1/2 



(32) 



Po, 



For an arbitrary k^ Hp we seek a solution of Equation (j30p in the form: 
Vz{t, X, z) = V{z) exp(At + ik±^ x). Introducing new variables: 



*(i?) == Vi?F(z), R{z) 



'^AO 



irns"p 



= ^/ffp 



(33) 
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we can rewrite Equation pop in the form of a 1-D Schrodinger equation for the 
function "^(R): 

where p = pp e^^'^'', wao = -Bo/\/Po i'' ^hc Alfven speed, the parameter a is 

and the potential U{R) has the fohowing asymptotic behavior: U{R — >■ 0) = 
k'^Hp/R^ and U{R — )■ oo) = a/R. For the existing of the instabihty, the potential 
U{R) should have a negative minimum. For example, for a long wavelength 
instability {k^Hp ^ 1) and when qpQ > 1, the potential U{R) has a negative 
minimum, and the instability can be excited. When the potential U{R) has a 
negative minimum and since U{R ^- 0) > and U{R —J- oo) > 0, there are two 
points i?i and R2 (the so-called turning points) in which U{R) = 0. Using the 
equations C/(i?i,2) = and Equations (|M)) - (|35[) we obtain the growth rate of 
the instability as 

_ /3. M„ns [RlR2{2 + Rl + R2)]'^^ , . 

H, (l + i?i)(l + i?2) ' ^-^^^ 

where we have used /3* = /3p-y/9pO- Note that Equation p6p is consistent with 
the simple estimate ([T2|). 
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